plot(-10:10, -10:10, xaxp=c(-10,10,20), yaxp=c(-10,10,20),type='n',xlab='',ylab='',asp=1)
xx <- c(-1,-1,0,0,1,2,3,3)
yy <- c(9,4,3,0,2,0,-2,-5)
cor(xx,yy)

points(xx, yy, col='blue', bg='blue', pch=21, cex=1.5)
abline(h=0)
abline(v=0)
grid()
beta <- (mean(xx*yy) - mean(xx)*mean(yy))/(mean(xx*xx) - (mean(xx))*(mean(xx)))
alpha <- mean(yy) - beta*mean(xx)

zz <- seq(-10,10,by=0.1)
bestfitline <- function(x) {
  alpha + beta*x
}
points(zz,bestfitline(zz),type='l',col='red', lwd=2)
